Frédéric’s mail

Salut Wang: tu peux mettre ce script ou l’adapter? Il leur faut aussi pouvoir telecharger le fichier de données (a mettre au meme endroit que le script!

J’imagine que la moitié ne sauront pas même faire ça… sans parler des problemes d’extension sous > windows; à prendre les gens pour des idiots on en fait des idiots)

The Plan of Work

  • Installing R: and package BN learn, test if the package is well
  • Dose-response modelling
    • KBrO3 -> GSH dataset
    • Préparer un script (html) :
    • 45 min de formation de R
      • Vecteur c()
      • Assignement x =
      • Concept
      • Fonction (demo)
      • Packages
    • Faire ajuster à la main un modèle de Hill avec quatre parametre à jouer Top Bottom

1 Preperation

1.1 Set Working Directory

R is always pointed at a directory on your computer.

  • R loads files (Plain text files / Data files etc.) from this directory.
  • R outputs files (Plain text files / Data files etc.) to this directory.

  • Set Working Directory (WD)
    • On Windows
      • in R interface choose : Menu - File - Change dir...
      • Choose a local directory (folder)
        • for example I use C:/Users/Admin/Documents/Konstanz_summer_school/Day2 as my working directory
    • On MacOs
      • in R insterface choose : menu - Misc - Change Working Directory...
        • for example I use ~/Documents/Konstanz_summer_school/Day2 as my working directory
    • Verify your WD
      • type getwd() in R console prompt: after >, then press enter key
      • R should print your current WD

1.2 Download dataset and load it to R

  • The dataset is accessible via the following two links Link 1 or Link 2
  • Click on one of two Links, a new web page will open then you have two options
    • Option 1: press Windows Shortcut : Control + S or MacOS Shortcut : Command ⌘ + S
    • Option 2: right click your mouse (or touchpad) and choose save as
  • save data_GSH.txt to your working directory (in my case /Documents/Konstanz_summer_school/Day2)
  • once this is done, type list.files() in R console prompt: after >, then press enter key
  • R should print all files in your working directory. Check if data_GSH.txt is listed
  • if not click on

Copy - Paste in R console prompt: after >, then press enter key

tmpData = readLines("https://antoinechn.github.io/Konstanz-SummerSchool/data/data_GSH.txt")
writeLines(tmpData, "data_GSH.txt")
rm(tmpData)
  • once this is done, type list.files() in R console prompt: after >, then press enter key
  • R should print all files in your working directory. Check if data_GSH.txt is listed
  • if still not
    • check your internet connection !!!
    • or contact us.
list.files()

1.3 Load data

This is GSH (after 1h)

data_GSH = read.table("data_GSH.txt", sep="\t", header=TRUE)

head(data_GSH)
#>   KBrO3_mM group pct_control
#> 1    0.000     1   91.521770
#> 2    0.375     1   63.591550
#> 3    0.750     1   50.508870
#> 4    1.500     1   14.285990
#> 5    3.000     1    2.477166
#> 6    6.000     1    0.377906
data_GSH = read.table("https://antoinechn.github.io/Konstanz-SummerSchool/data/data_GSH.txt", sep="\t", header=TRUE)
# data_GSH = read.table("https://trello-attachments.s3.amazonaws.com/59552633a56853b5f6d1f1f9/5a92760183331b90626ce166/678efebfc21386538df65c23120fc741/data_GSH.txt", sep="\t", header=TRUE)

head(data_GSH)

2 Data visualisation using plot() function

Plot the data to get a sense of it

plot(data_GSH$KBrO3_mM, data_GSH$pct_control)


You have plenty of options, the most useful are:

# nicer
plot(data_GSH$KBrO3_mM, data_GSH$pct_control, las=1)
# professional
plot(data_GSH$KBrO3_mM, data_GSH$pct_control, las=1,
     xlab="KBrO3 (mM)", ylab="GSH decrease (%)")
# excellent
plot(data_GSH$KBrO3_mM, data_GSH$pct_control, las=1,
     xlab="KBrO3 (mM)", ylab="GSH decrease (%)", main="My Beautiful Data")
# clever
plot(data_GSH$KBrO3_mM, data_GSH$pct_control, las=1,
     xlab="KBrO3 (mM)", ylab="GSH decrease (%)", main="My Beautiful Data",
     log="y") # you can use log="x" or "y" or "xy" or "" 
# genius level
plot(data_GSH$KBrO3_mM, data_GSH$pct_control, las=1,
     xlab="KBrO3 (mM)", ylab="GSH decrease (%)", main="My Beautiful Data",
     log="y", pch=10, cex= 1.3, col="red") # set the type of mark, size, color

# etc. see help(par) or ?par for all the parameters you can set
# nicer
plot(data_GSH$KBrO3_mM, data_GSH$pct_control, las=1)

# professional
plot(data_GSH$KBrO3_mM, data_GSH$pct_control, las=1,
     xlab="KBrO3 (mM)", ylab="GSH decrease (%)")

# excellent
plot(data_GSH$KBrO3_mM, data_GSH$pct_control, las=1,
     xlab="KBrO3 (mM)", ylab="GSH decrease (%)", main="My Beautiful Data")

# clever
plot(data_GSH$KBrO3_mM, data_GSH$pct_control, las=1,
     xlab="KBrO3 (mM)", ylab="GSH decrease (%)", main="My Beautiful Data",
     log="y") # you can use log="x" or "y" or "xy" or "" 

# genius level
plot(data_GSH$KBrO3_mM, data_GSH$pct_control, las=1,
     xlab="KBrO3 (mM)", ylab="GSH decrease (%)", main="My Beautiful Data",
     log="y", pch=10, cex= 1.3, col="red") # set the type of mark, size, color



3 Modelling

3.1 Hill dose-response model

Let’s get serious:

  • define the Hill dose-response model as a function
  • let’s try it
  • plot it…
# define the Hill dose-response model as a function
Hill = function (x, EC50, n) {
  return(x^n / (x^n + EC50^n))
}

# let's try it
Hill(x=10, EC50=1, n=1)

# plot it...
x = seq(0, 10, 0.1)
y = Hill(x=x, EC50=1, n=3)
plot(x, y, type="b")
# define the Hill dose-response model as a function
Hill = function (x, EC50, n) {
  return(x^n / (x^n + EC50^n))
}
# let's try it
Hill(x=10, EC50=1, n=1)
#> [1] 0.9090909
# plot it...
x = seq(0, 10, 0.1)
y = Hill(x=x, EC50=1, n=3)
plot(x, y, type="b")

Problem:

  • it goes up from zero to one
  • and we want it to go down from 100 to 0

3.2 Extended Hill dose-response model

So we change it a bit

  • define an extended Hill dose-response model as a function
  • try it…
  • plot it…
# define an extended Hill dose-response model as a function
Hill = function (x, EC50, n, from, diff) {
  return(from + diff * (x^n / (x^n + EC50^n)))
}

# try it
Hill(x=10, EC50=1, n=1, from=1, diff=10)

# plot it...
x = seq(0, 10, 0.1)
y = Hill(x=x, EC50=1, n=3, from=1, diff=10)
plot(x, y, type="b")
# define an extended Hill dose-response model as a function
Hill = function (x, EC50, n, from, diff) {
  return(from + diff * (x^n / (x^n + EC50^n)))
}

# try it
Hill(x=10, EC50=1, n=1, from=1, diff=10)
#> [1] 10.09091
# plot it...
x = seq(0, 10, 0.1)
y = Hill(x=x, EC50=1, n=3, from=1, diff=10)
plot(x, y, type="b")


Question

  • Can you find values of from and diff that would match the data range?
    • first: what is the data range?
    • assign the values found to variables A and B
# assign the values found to variables A and B
A =
B =

now replot the data and overlay with the model

Note : make sure to define the axes limits before doing overlays…

xlims = c(0, 6) # a vector of min and max
ylims = c(0, 120)
plot(data_GSH$KBrO3_mM, data_GSH$pct_control, xlim=xlims, ylim=ylims, col="red")
par(new=T) # will overlay the next plot
x = seq(0, 6, 0.1)
y = Hill(x=x, EC50=1, n=3, from=A, diff=B)
plot(x, y, type="l", xlim=xlims, ylim=ylims, xlab="", ylab="") # turn off labels

now play with EC50 and N to fit the data better

  • A = 100
  • B = -100
  • EC50 = 0.5
  • N = 1.8

3.3 Notes

There are many tools in R that can be used to fit automatically your model one of them is PROAST (from the RIVM) but we can also access PROAST online (more fun)